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1 Introduction 

During the past decades, there has been an increasing fascination and sur- 
prises with diverse quantum many-body effects. With the magical touch of in- 
teraction a simple electron system may assume insulating, metallic, magnetic 
or superconducting states according as the control parameters are changed. 
Strongly correlated electron systems, as exemplified by the high-Tc supercon- 
ductors and their host materials realized in transition-metal oxides, as well as 
by organic metals, have provided us with an ideal playground, where various 
crystal structures with band-filling control and band- width control etc provide 
the richness in the phase diagram [39] . 

On the other hand, there is a long history of the interests in non- equilibrium 
phase transitions. Statistical mechanically, there is an intriguing problem of 
how we can generally define the notion of a "phase" in non-equilibrium sys- 
tems, but we can still discuss individual systems in specified non-equilibrium 
conditions to extract more general viewpoints. Now, if we combine the above 
two ingredients, namely, if we consider strongly- correlated electron systems in 
non- equilibrium, we plunge into an even more fascinating physics. In fact, re- 
cent years have witnessed an upsurge of interests in non-equilibrium states in 
many-body systems with drastic changes in the electronic states in strong dc 
electric fields, in intense laser fields, etc. 

Developments in fabrication techniques such as realization of clean thin 
films with electrodes attached have triggered several groundbreaking experi- 
ments, e.g., non-linear transport measurements in thin films [1, 60, 19, 50, 2, 4], 
in layered systems [24 and observations of clean metallic states in heterostruc- 
tures [3|. Non- linear phenomena in correlated electron systems now begin to 
attract interests in a wide range of researchers: One obvious area of appli- 
cation is future-generation electronic devices, where a high sensitivity of a 
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Fig. 1. Non- linear transport and optical response 



system near a phase boundary to external conditions may lead to drastic 
functionalities [I]. However, even more attractive is its relevance to fundamen- 
tal physics, especially, to non-equilibrium statistical physics, where we can 
observe the behavior of various phase transitions taking place under non- 
equilibrium conditions. 

The purpose of the present article is to discuss the nonequilibrium metal- 
insulator transition in strongly correlated electron systems [57, 58, 7, 8 , which 
is known, for equilibrium systems, as Mott's transition. Before going into de- 
tail, we first give a brief introduction of the transition, and discuss how quan- 
tum breakdown through non-adiabatic transitions in nonequilibrium becomes 
relevant in non-linear transports. 

The model we study is the single-band Hubbard model, which is the sim- 
plest possible one that captures many essential properties of correlated elec- 
tron physics. The Hamiltonian reads 



where Ci a annihilates an electron on site i with spin a, rii a = c ia .Ci a the 
number operator, U the strength of the on-site Coulomb repulsion, and fhop 



the hopping integral. The filling n = — 2_,( n i] with L the number 



of sites, is an important parameter, which changes the groundstate property 
drastically. 

When the band is half-filled with one electron per site on average (n = 1), 
each electron tends to be localized on a separate lattice site and the spin 
tends to be antiferromagnetically correlated. When U/t is large enough, the 
groundstate is insulating, which is called Mott's insulator (Fig[2|a)), and the 
groundstate is separated from the charge excited states with a many-body 
energy gap — Mott gap. When we inject carriers (usually with a chemi- 
cal doping by adding or replacing to other elements) to increase (electron 
doping) or decrease (hole doping) the filling from unity, the Mott gap col- 
lapses for large enough doping, and the system becomes metallic. This is the 
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Fig. 2. Metal-insulator transition in equilibrium due to doping: (a) A Mott insulator 
realized at half filling, (b) A hole-doped metallic state, (c) An electron-doped metallic 
state. 



metal-insulator transition or the Mott transition, which is widely observed in 
strongly correlated materials. In these materials, a state occupied simultane- 
ously by up an down-spin electrons — which we call a doublon — and holes 
carry the current. After the discovery of the high-temperature superconduc- 
tivity in cuprates, carrier-doped Mott insulators have been subject of a huge 
number of experimental and theoretical studies. 
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Fig. 3. (a) Schematic experimental configuration, (b) Carriers (doublons and holes) 
created by an external electric field. 



Now, let us consider what will happen if we attach a set of electrodes 
to a strongly-correlated sample, and apply a large bias voltage across the 
electrodes (Fig. [3] (a)). Although the setup may seem simple enough, there is 
a profound physics involved. 

In regions near the electrodes (or near the interface in the case of het- 
erojunctions between a strongly-correlated and ordinary materials), a "band 
bending" similar to doped semiconductors can take place and lead to an in- 
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terface Mott transition when the filling becomes one [55] • The width of the 
insulating layer changes as the applied bias is changed, which dominates the 
behavior of the non-linear transport (i.e., the I—V characteristics). The result 
(with DMRG + Hartree potential) for the band-bending effects in this case 
can be understood if we assume a local equilibrium for the relation between 
the density of electrons and the potential. The local properties are determined 
by the Hartree potential governed by Poisson's equation, which in turn deter- 
mines the local chemical potential and controls the metal-insulator transition. 

Even more interesting, however, is the case where we no longer have local 
equilibrium. Specifically, a quantum many-body breakdown of a Mott's insu- 
lator takes place when the applied electric field is large enough and creates 
doublons and holes in the Mott insulating groundstate (Fig.[3](b)) |58, ]Hl [57] . 
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Fig. 4. Dielectric breakdown of a Mott insulator in a strong electric field due 
to many-body Landau-Zener transition: The groundstate and excited states with 
charge excitations are separated by an energy barrier, and quantum tunneling among 
many-body states takes place when the electric field is strong enough. 



The creation mechanism is a many-body analog of the "Zener breakdown" , 
well-known in semiconductor physics |68| . Namely, while we cannot use the 
notion of the electronic band structure for correlated electron systems, we can 
envisage the carrier-creation process as a tunneling across a kind of barrier. As 
displayed in Fig. [4l production of carriers occurs through tunneling between 
the Mott insulating groundstate and excited states with doublons and holes. 
If we denote the distance between a doublon and a hole by l^h, the energy 
profile as a function of Idh roughly reads 
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AE~U- l dh F, (2) 

where F is the strength of the electric field. The profile curve reaches the 
energy before the creation of the doublon-hole pair for the separation at which 

AE~U- l dh F = 0, (3) 

to which the tunneling becomes possible. There is a threshold field strength for 
this process to occur. This is because larger quantum fluctuations are required 
to have a larger separation between doublon and hole in the Mott insulator. 

In other words, the overlaps of many-body wave function of the ground- 
state and excited states decrease rapidly for large l d h- We can formulate this 
with the Landau-Zener picture in the time-dependent gauge for the external 
electric field, for which, as we shall show below[57](eq. ([5T1l ). the threshold field 
strength is given by 

Ftt = f£, (4) 

where A C (U) is the charge gap (i.e., the Mott gap), and the tunneling rate 
per length is given by 

IF 

r(F)/L = -—aln[l-p(F)], (5) 

where p(F) — e~~ 7T ^ r ~ is the tunneling probability and a a non-universal con- 
stant depending on the detail of the system. The tunneling rate r(F)/L, 
being related to the production rate of carriers, is directly related to physical 
properties in the bulk if the interface effect is neglected. Indeed, such a non- 
linearity in the I — V characteristics has been observed in real materials, most 
prominently in a one-dimensional copper oxide [41] . 

Another interesting consequence of eq.((3|) is that it gives the "critical sep- 
aration" of the doublon-hole excitation l d h = U/Fth- This has to do with the 
convex shape of the energy profile against Idh (FigH]), which is reminiscent 
of the energy profile for the standard nucleation theory that treats the crit- 
ical size of a stable-phase droplet to grow without being crushed, although 
the physics involved is quite different. In the present case, when the field 
is greater than the threshold, the electric-field induced metallic state, where 
doublon-hole pairs continue to be created, becomes the stable phase. 

The first goal of this article is to derive the relations presented above and 
study the creation mechanism of carriers (this part is the extended argument 
of our papers [581 157j ) . We need to treat the process quantum mechanically 
and in a many-body formulation. In doing so, we present a renewed and unified 
interpretation of the Zener transition of insulators. The key quantity is the ef- 
fective Lagrangian of quantum dynamics (see H2.2l for a detailed introduction) 
which is define by [57] 

£(F)=-^ t limilnE(t), (6) 
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where S(t) is the groundstate-to-groundstate transition amplitude and L d is 
the volume of the d dimensional system and L the linear size. There is a deep 
relation between the theories of dielectric breakdown in condensed matter and 
a branch in quantum field theory known as non-linear quantum electrodynam- 
ics (QED) (table: [J). The effective Lagrangian defined above coincides with 
the Heisenberg-Euler effective Lagrangian for non-adiabatic evolution [20]I57|. 
The effective Lagrangian have been used to study the Schwinger mechanism 
of electron-positron pair production from the QED vacuum in strong electric 
fields 52J. In fact, we show that the Schwinger mechanism and the Zener 
tunneling are equivalent, where the effective action coincides if we consider 
the breakdown of simple Dirac type band insulators. Furthermore, the effec- 
tive action gives the non-adiabatic extension of the Berry phase theory of 
polarization. 





dielectric breakdown in cond. matter 


non-linear QED 


mechanism 

excitation 

effect of interaction 

non-linear polarization 


Zener breakdown [58] 
electron (doublon)-hole pair 
many-body Landau-Zener 
cross correlation (ME effect) 
non-adiabatic Berry phase theory [57] 


Schwinger mechanism [52] 
electron-positron pair 
back reaction 
photon-photon interaction 



Table 1. Relation between the theory of dielectric breakdown in condensed matter 
and non-linear QED from the point of view of the effective Lagrangian. 



In the latter part of the article we shall discuss the effect of annihilation 
of doublon-hole pairs (Fig. [5|. In a one-dimensional system, a doublon and 
a hole cannot pass each other without being pair-annihilated even as virtual 
processes. Since the groundstate is locally stable, the many-body state tends 
to remain in the ground state, but there should be a finite probability for 
the state to "branch into" excited states through many paths in the many- 
body energy space. Thus, the long-time behavior of the wave function involves 
numerous scattering processes in the energy space, where, as we shall see, the 
phase interference plays a key role. We can indeed regard the phase before the 
dielectric breakdown takes place as a dynamical localization in the many-body 
energy space, which reduces the tunneling rate and makes the groundstate 
survive [8] [7]. A statistical mechanical treatment helps in understanding this, 
and we briefly discuss it in terms of the quantum walk. 

A brief comment on the numerical methods used in the article: In order 
to understand the non-equilibrium processes, we need to integrate the time- 
dependent, many-body Schrodinger equation to look at the evolution of the 
many-body wave function in, say, the Hubbard model in strong electric fields. 
This is a formidable task, for which no analytically exact treatment is known, 
so that we rely on several numerical methods, which include the exact di- 
agonalization and the time-dependent density matrix renormalization group 



Title Suppressed Due to Excessive Length 7 



<$> O <$><$> ft <$> Electric field F 




separation 

Fig. 5. Annihilation processes for carriers in a correlated electron system. 

method (td-DMRG; The version of td-DMRG we adopt is the one proposed 
by White and Feiguin[63]). 

2 Non-adiabatic evolution and pair creation of carriers 
2.1 Electric fields and gauge transformation 

When we describe a system in finite electric fields, we can choose from two 
gauges. One is the case where we have a slanted electrostatic potential, with 
the gauge field = (Fx,0) for a one-dimensional system with F = eE 
being the electric field, while the other represents the electric field via a time- 
dependent vector potential, = (0, —Ft). 




Fig. 6. (a) Time-independent gauge, (b) Time-dependent gauge. 



In the first gauge, the tilted potential enters in the Hamiltonian as 

H(F) = H Q + FX, X = ^jn„ (7) 

3 
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where X is the position operator [JSj and H the original Hamiltonian. This 
gauge is incompatible with systems with periodic boundary conditions. The 
Hamiltonian H(F) in fact becomes an unbounded operator in an infinite sys- 
tem, since one can lower the energy indefinitely by moving an electron to 
j -> -oo. 

In the other gauge, which we call the time-dependent gauge, the hopping 
term of the Hamiltonian becomes time-dependent as 

H (#)) = -t hop (^ (t) + h-c) + V, (8) 

where <P(t) represents a time-dependent Aharonov-Bohm(AB) flux, 

4>{i) = $(t)/$ = FLt/h. (9) 

Physically, this gauge amounts to considering a periodic system (a ring) and a 
magnetic flux piercing the ring, where the time-dependent flux induces electric 
fields by Faraday's law (Fig. [6]). For a higher-dimensional system, the ring 
becomes a (generalized) torus. The time-dependent gauge is suited for periodic 
systems since it is compatible with the lattice translation symmetry. The 
electric current operator is obtained by differentiating the Hamiltonian by A 1 
as 

JW = = _i * hop ? (? 3t *^«° - ^^cU a c l+la ) . (10) 



There exists an important operator relation among H, J and X, 

nm,x], (ii) 



which comes from Heisenberg's equation of motion for the current operator, 

J(0) = j- t x. 

We can relate the two gauges with a twist operator [S] defined by 

S(0) = e-**** (12) 

and the two Hamiltonians are related by a gauge transformation generated by 
the twist operator, i.e., 

H(F) - gi(d>(t))H(4,(t))g(<t>(t)) ~ ^ t (0W)^. 9 (0W). (13) 



2.2 Heisenberg-Euler effective Lagrangian 

We first discuss the non-adiabatic evolution of electron wave functions in insu- 
lators (either one-body or many-body) in strong electric fields. Let us consider 
an insulator at T — and F = 0, which is described by the groundstate wave 
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function \&o). We then switch on the electric fields at t = to study the 
quantum mechanical evolution of the system. We limit our discussions to co- 
herent dynamics and ignore the effect of dissipation due to heat bath degrees 
of freedom as well as boundary effects near the electrodes. 

A key quantity to study the non-adiabatic evolution and quantum tun- 
neling in strong electric fields is the (condensed-matter counterpart to the) 
effective Lagrangian introduced for QED by Heisenberg and Euler[2Q]. In the 
time-independent gauge, the electrons are described by the solution of the 
Schrodinger equation, 

\#(t))=e- aH W\9 ), (14) 

where we have put h = 1. The overlap of the solution with the groundstate for 
F = — groundstate-to-groundstate transition amplitude — should contain 
the information on the tunneling processes, so we define 

S (t) = {V \e- UH{F) \Vo)e ltEo , (15) 

where we have factored out the trivial dynamical phase of the groundstate, 
E = (<P\H(F = 0)\\P). In the case of the time-dependent gauge, we need 
to be careful, since the groundstate is </> dependent. If we denote \0;<fi) as the 
instantaneous groundstate of H(<fr), the groundstate-to-groundstate transition 
amplitude becomes 

S(t) = (O;0(r)|f e "^/o Tff ^ s » ds |O;0(O))e^r £o ^» ds , (16) 

where T stands for the time ordering, and Eo((f>) — (0; (f>(r)\H (0)|O; </>) the 
dynamical phase of the groundstate. 




Fig. 7. The original problem studied by Callan and Coleman in which quantum 
tunneling from an unstable vacuum is considered [23] . 



We define [57] the effective Lagrangian by 

£(F ) =-^ t HmilnE(t ) , 



(17) 
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where L d is the volume of the d dimensional system with a linear dimension 
of L. We can also regard the Lagrangian as the exponent of the asymptotic 
behavior of the amplitude, S{r) ~ e lrL £ ( F ). When the electric field is large 
enough, the groundstate becomes unstable with the quantum tunneling to 
excited states activated. The tunneling rate is described by the imaginary 
part of the effective Lagrangian, 

r(F)/L d = 2lmC{F), (18) 

which gives the rate of the exponential decay of the vacuum (groundstate). 
In the quantum field theory, the decay rate of an unstable vacuum has been 
discussed by Callan and Coleman, where the tunneling takes place when the 
potential is suddenly changed by an external field 23J (Fig. [7]). 




real space 

Fig. 8. Two models of the dielectric breakdown studied by Zener. (a) Time- 
independent gauge, (b) Time-dependent gauge. 

As we shall see later in several models, in the theory of dielectric break- 
down, the tunneling corresponds to creation of charge carriers. In band insu- 
lators the carriers are electrons and holes, while in Mott insulators they are 
doublons and holes. If we neglect boundary effects and assume that all the 
carriers are absorbed by electrodes, we can conclude that the tunneling rate 
is proportional to the leakage current, i.e., 

J lcak cx r(F)/L d . (19) 

Indeed, this is the original picture of Zener when he calculated the leakage 
current in a simple band insulator [8] Zener has studied the dielectric break- 
down in a simple one-dimensional insulator using the time- independent gauge 
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[59"] as well as the time-dependent gauge [55]. In the former, he has calcu- 
lated the tunneling probability of Bloch functions in constant electric fields 
to obtain the tunneling rate. In the time-dependent gauge, he has described 
the problem as a system with a time-dependent Hamiltonian represented by a 
two-by-two matrix to study the tunneling near the level anti-crossing, which 
is now known as the Landau-Zener transition [38] I68j. 

The reason why we have called C(F) the effective Lagrangian is that it 
coincides with the Heisenberg-Euler effective Lagrangian in QED [20 . They 
have studied the dynamics and nonlinear responses of the QED vacuum in 
strong electric fields by calculating the effective Lagrangian (for a review, 
see e.g. [15]). By integrating out high-energy degrees of freedom (polarization 
processes due to electron-positron creation/annihilation) they arrived at an 
effective description of the low-energy degrees of freedom, namely the quan- 
tum correction, originating from the fluctuation of the QED vacuum, to the 
Maxwell theory of electromagnetism. Indeed, if we apply our formalism to 
band insulators with Dirac-type (mass-gapped) dispersions, the effective La- 
grangian coincides with the Heisenberg-Euler Lagrangian with some modifi- 
cations coming from the Brillouin zone structure of the Bloch waves as will be 
shown in the next section. The correspondence between the two phenomena is 
straightforward: The ground-state of the insulator translates to the QED vac- 
uum, charge excitations to the electron-positron pairs. The tunneling rate also 
has its QED counterpart, namely the vacuum decay rate due to the Schwinger 
mechanism — creation of electron-positron pairs in strong electric fields 52J. 

Related theories 

We conclude this section with comments on the relation of the effective La- 
grangian approach to earlier theoretical frameworks. 

Berry's phase theory of polarization: 

In the Berry's phase theory of polarization [47] [29] [48] [49] [44] , the ground- 
state expectation value of the twist operator e~ %2 ^ x , which shifts the phase 
of electron wave functions on site j by — =j-j [44] , plays a crucial role. It was 
revealed that the real part of a quantity 



gives the linear-response electric polarization, P e \ = — Kew [48], while its imag- 
inary part gives a criterion for metal- insulator transition, i.e., D — 47rlmw; is 
finite in insulators and divergent in metals [49] . The present effective action 
is regarded as a non-adiabatic (finite electric field) extension of w. To give a 
more accurate argument, recall that the effective Lagrangian can be expressed 
as 



w = 



-^l n <0|e-^|0) 

Z7T 



(20) 




(21) 
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for d = 1. Let us set r = h/LF and consider the small F 1 limit. For insulators 
we can replace H with the groundstate energy E to have C(F) ~ W F in the 
linear-response regime. Thus the real part of Heisenberg-Euler's expression [2U] 
for the non-linear polarization Pre(F) = —dC(F)/dF naturally reduces to 
the Berry's phase formula P e \ in the F — > limit (cf. ea.(|27[) below). Its imag- 
inary part, which is related to the decay rate as ImPHE(F) = —^ ^qJ? , 
reduces to —D/An and gives the criterion for the transition, originally pro- 
posed for the zero field case. 

Non-Hermitian quantum mechanics 

The dielectric breakdown of Mott insulators was also studied in the framework 
of non-Hermitian quantum mechanics [1 7| [67] ■ Fukui and Kawakami studied a 
non-Hermitian Hubbard model in which the leftward and rightward hopping 
integral are assumed to be unequal [17] . The non-Hermiticity is assumed to 
represent the coupling of the system with a "dissipative environment" . With 
the Bethe ansatz solution they have observed the gap between the groundstate 
and the first excited state to close when the hopping asymmetry is large 
enough. It seems that the remaining question is to relate this result with 
measurable quantities. 

2.3 Zener breakdown of band insulators revisited — 
Non-adiabatic geometric phase and the Schwinger mechanism 





(a) (b) 

Fig. 9. (a) Energy levels of a Dirac band in 2D. (b) A one dimensional slice of the 
higher dimensional Dirac band in which carriers (doublons and holes) are created 
by an external electric field in the ku direction. 
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Before examining the dielectric breakdown of correlated electron systems, 
let us first discuss the dielectric breakdown of band insulators in an electric 
field F within the effective-mass picture. This will turn out to be heuristic, 
since we can obtain an analytic expression for the effective Lagrangian which 
can be readily applied to general band insulators. 

For simplicity we take a pair of hyperbolic bands £±(k) = zty/V 2 + v 2 k 2 
(considered here in d spatial dimensions), where 2V is the band gap, — (+) 
denote the valence (conduction) band, and v the asymptotic slope of the 
dispersion. 

We first obtain the groundstate-to-groundstate transition amplitude with 
the time-dependent gauge in the periodic boundary condition. There, a time- 
dependent AB-flux in units of the flux quantum, 4>{t) = FLr/h (with the 
electronic charge e = 1 and L being the system size), is introduced to induce 
an electric field F, which makes the Hamiltonian time dependent as 

H(Hr))= J2 £a(k+|V)e||) C t(k) CQ (k). (22) 

k,a=± ^ ' 

Here eii is the unit vector parallel to F, and cjj(k) the creation operator with 
spin indices dropped. If we denote the ground state of -ff (</>) as |0; </>) and its 
energy as Eo(4>), the groundstate-to-groundstate transition amplitude reads 

3{t) = (O;</)(T)|fe-^/o Tff W s » ds |O;0(O))e^r £o ^» ds , (23) 

where T stands for the time ordering. The effective Lagrangian C(F) for 
the quantum dynamics is defined from the asymptotic behavior, 5(r) ~ 

e %TL d C{F)_ 

The dynamics of the one-body model can be solved analytically (Fig[9] 
(b)), since we can cut the dispersion in d spatial dimensions into slices, each 
of which reduces to Landau-Zener's two band model in lD[38l |68]. Namely, 
if we decompose the k vector as (kj_,fc||), where k^ (fc||) is the component 
perpendicular (parallel) to F, where each slice for a given k^ is a copy of 
Landau-Zener's model with a gap Z\b an d(k) = 2y/V 2 + v 2 k\. The Landau- 
Zener transition takes place around the level anti-crossing for which k\\ + 
X"</>(t) moves across the Brillouin zone(BZ) in a time interval St = h/F. 
The process can be expressed as a scattering and the Bogolubov coefficients 
between the "in" and "out" states (see Figj9]) is given by the solution to the 
two band problem, i.e., 

cj.(k) - Vl-p(k)e-*< k ><4(k) + Vp(k)c f _(k), 

cl(k) -» -Vp(k)4(k) + v/l-p(k)e ix ( k )cL(k). (24) 

Here the tunneling probability for each k is given by the Landau-Zener(LZ) 
formula[3EllE], 

(Z\ band (k)/2) 2 1 



p(k) = exp 
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On the other hand, the phase x(k) = — #(k) + 7(k) appearing in the Bogol- 
ubov coefficients consists of the trivial dynamical phase, H9(k) = J Q T e+(k + 
^<j>{s)e\\)ds, and the Stokes phase [551 125] . 



7 (k) = -Im 



G?S- 



,-i(Zl band (k)/2) 2 S 



cot(v.Fs) 



1 



(26) 



The Stokes phase, a non-adiabatic extension of Berry's geometric phase [10] . 
depends not only on the topology of the path but, unlike the adiabatic coun- 
terpart, also on the field strength F [26]. In terms of the fermion operators 
the groundstate is obtained by filling the lower band |O;0) = rik c -(k~ 
^S(/>e||)|vac), where |vac) is the fermion vacuum with c±(k)|vac) = 0. If we as- 
sume that excited charges are absorbed by electrodes we obtain from eqs. ([2"5]l . 
121 



-F 



Re C{F) 

Im C(F) =-F [ 
Jb 



dk 7 (k) 



bz (27r) d 2ir 
dk 1 



In[l-p(k)], 



/ BZ (2ir)d Air 

where the dynamical phase 8 cancels the factor e R Jo E o(<f>( s )) ds m e g (J23 



(27) 




F/R 



th 

Fig. 10. The dependence of the conductivity on the electric field in the non-linear 
regime for band insulators with spatial dimension d — 1,2, 3. The inset zooms in the 
threshold region. 



Integration over k in eg. ([27]) leads to the groundstate decay rate per vol- 
ume for a d-dimensional hyperbolic band, 
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F ( F 



) 



(d-l)/2 



r(F)/L d 



{2TT) d ^h \ V 





-, d-1 



(d+l)/2 e 




erf 



(28) 



The threshold for the tunneling is governed by the most nonlinear (actually 

essentially singular) factor in the above formula, namely e~^ n ^ F , so that the 
threshold electric field is given by 



Although an analytic integration eq. (|27p is possible for a Dirac band (= hyper- 
bolic valence and conduction bands) , the expression is valid for general band 
dispersions. In fig-EH the leakage current divided by the field strength, which 
is proportional to r(F)/FL d , is plotted for the spatial dimension d = 1, 2, 3. 
The independence is essentially 



which has a threshold behavior as shown in the inset of the figure. Above the 
threshold, two regimes exist. In the medium field regime, the current scales 
as Jieakage ~ F^ d+1 ^ 2 where the power depends on d. However, when the field 
strength is even stronger, the error function appearing in eq. (|28p . which is 
due to the lattice structure (with the k integral restricted to the BZ), starts 
to take an asymptote (erf(x) ~ (2/y / 7r)x). Then various factors (including a 
power of F) cancel with each other, and the leakage current in the F — > oo 
limit approaches a universal function, 



where the d dependence disappears up to a trivial (/-dependent numerical 
factor. This prediction on the non-linear transport can be checked experimen- 
tally including low-dimensional systems such as carbon nanotubes (d = 1). 
Graphene (d = 2) is also interesting, but this system has a massless Dirac 
dispersion, so that a special treatment should be required. 

Comparison to Heisenberg-Euler-Schwinger's results in non-linear QED 

Let us have a closer look at the decay of the QHE vacuum. In 1936, Heisen- 
berg and Euler studied Dirac particles in strong electric fields, and discussed 
non-linear optical responses of the QED vacuum — vacuum polarization — 
in terms of an effective Lagrangian [20] . Later, Schwinger refined their ap- 
proach and calculated the vacuum decay rate [52^1 . Up to the one-loop level, 

1 For references on the effective- act ion approach of non-linear electrodynamics, see 




(29) 



^leakage OC F ^ 6 F , 



(30) 




(31) 



[nun]- 



16 Takashi Oka, Hideo Aoki 



Schwinger calculated the vacuum-to-vacuum transition amplitude using the 
proper time regularization method to obtain 



1 f°° ds_ 

2 



Fcot(Fs) - - 
s 



(32) 



for (3+l)-dimensional QED, where m e is the electron mass. The integrand 
has a pole in the complex domain and has an imaginary part, which gives 

™~/* -£§;?■*(- W> 1331 

the famous Schwinger's formula for the electron-positron pair creation rate 
[52] . where a = 1/137 is the fine-structure constant. 

Thus the expression for the QED effective Lagrangian, eq. ([32l) . coincides 
with the Stokes phase for the non-adiabatic Landau-Zener tunneling, except 
for a difference in the momentum integral. As we have mentioned above, an 
important difference in lattice systems is that the momentum integral is lim- 
ited to the Brilliouin zone, and the decay-rate acquires an extra factor (com- 
pare eq. pS)) with erf with eq. ([3"3"l) ). This modification changes the strong field 
limit of the leakage current which leads to the universal expression feq. (l3T10 . 
Another important difference, which is quantitative, appears in the thresh- 
old voltage: The threshold for band insulators -E t b h and = F^ nd /e = V 2 /vae 
(a: lattice constant) is many orders smaller than the threshold for the QED 

2 3 

instability E® EB = ^f - ~ 10 16 V/cm. For example, if we have an insulator 
with parameters a — 10~ 7 cm, v — 2thop = leV, V — leV, then the threshold 
becomes as small as E^ nd = 10 7 V/cm. 

Heisenberg and Euler's original aim was to discuss non-linear optical prop- 
erties of the vacuum in terms of AC. In fact, they calculated the effective 
Lagrangian in the presence of both electric and magnetic fields [20] , and ob- 
tained 

AC QEB {F) = C E2 - B2 + J^L [(E 2 - B 2 ) 2 + 7(E • B) 2 1 + . . . , (34) 
2 45m| L J 

where C is a diverging constant that we drop after renormalization. The elec- 
tric polarization can be obtained from the real part of the effective action 
via 

AP{F) = ~§F AC{F) - (35) 



If we plug in eq (|34[) , the non-linear polarization of Dirac particles becomes 

AP = ^—r (-AB 2 E + UB?,E + AE 3 ) + (36) 
45m* V 11 / 

oo 

= J2 p(n) (B)E n , (37) 



71=1 
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where Bu is the component of B parallel to E, and P"(B) the n— th or- 
der non-linear polarization. Thus we can examine nonlinear polarizations and 
cross correlations (a combined effect of E, B) with the effective Lagrangian, 
as touched upon in Table [1] 

2.4 Dielectric breakdown in a Mott insulator — many-body 
Landau-Zener transition and a nonequilibrium phase diagram 




0.0 0.2 0.4 0.6 0.8 1.0 0.0 0.2 0.4 0.6 0.8 1.0 



Mott Insulator 

U=1.5, V=0.0 




0.0 0.2 0.4 0.6 0.8 1.0 



Fig. 11. Adiabatic spectrum E„(4>) for a finite system (L — 10 here) obtained 
by the Lanczos method. We plot low-lying excitations in the half-filled subspace 
Nf = Ni = L/2. (a) a noninteracting system in a free space, (b) a band insulator 
(U/t hop = 0, V/t hop = 0.3), and (c) a Mott insulator (U/t hop = 1.5, V/thop = 0). 
The circles indicate avoided level crossings. 



Before applying the effective Lagrangian approach to the dielectric break- 
down of Mott insulators, we need to examine the excitation spectra, which is 
displayed in fig[TT] There we plot, for the half- filling, the many-body energy 
levels of the Hamiltonian, 

hw)) = -wE ( e<¥ * (t) 4u«rCfcr +h-c.) + uj2^i + vj2(-iym 

i(7 i i 

Here U is the Hubbard repulsion, V a staggered potential to introduce valence 
and conduction bands, so that U = V = corresponds to a noninteracting 
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system in a free space, U = 0, F / a band insulator, and a large U and 
V = a Mott insulator. In the figure, we have only plotted charge excitations 
(where the charge rapidities are excited in the language of Bethe-ansatz solu- 
tion) . As can been seen, levels cross in the free model while in the band and the 
Mott insulators an energy gap separates the ground state from excited states. 
The gap is 2V for the band insulator. The Hubbard Hamiltonian eq. (|38|) with 
U 7^ 0, V = is also exactly solvable in ID. Woynarovich used the Bethe 
ansatz method [65l [66] to study the ground state as well as the excited states 
(see also [551 GUI The two solid lines in Mott insulator's spectrum cor- 
respond to the ground state and the charge-excited state with one pair of 
complex charge rapidities, quantum numbers appearing in the Bethe-ansatz 
solution. The energy gap AE(U) between these states are known to converge 
to the Mott gap Am u(U) in the limit of infinite system, 

AE{U) -> A C (U). (39) 

An important feature in the spectrum of the Mott insulator is that level 
repulsion occurs at many places over the excited states. The repulsion comes 
from Umklapp electron-electron scattering, i.e., a scattering process in which 
the momentum sum changes by reciprocal lattice vectors. In band insulators 
level repulsions obviously come from one-body scattering as we have seen 
above. 

Why have we first looked at the adiabatic spectrum? There is an important 

relation between the adiabatic energy and the current expectation value. From 

dE (WldH/dXlW) 

the Hellmann-Feynman theorem, i.e., — = ; — ■ — ; for H(X)\W(X)) = 

dX (&\&) Wl K " 

E{X)\W(X)), we obtain 

J n {<t>) = {n-^\J^)\n-4>) 

which is valid for all <f>. If we expand it around (f> — 0, we get 

L \ d 2 E n (0) 



U<t>) - MO) + [-1 + OWO. (41) 

Using (f> = FLt/h and defining the transport coefficients T> n by J n ((f>) = 
J n (0) + V n Ft + 0{F 2 ), we obtain 

When we focus on a dissipationless adiabatic transport at T = 0, the total 
current thus reads 

(J(i)) =V Q (L)Ft, (43) 



Title Suppressed Due to Excessive Length 



19 



which is determined by the Drude weight (charge stiffness) T> a (L). As we can 
see in FigfTTl even for insulators ((b) and (c)), the Drude weight T> (L) of a 
finite system is not necessarily zero. If we remember Kohn's criterion\61\ for 
metal-insulator transitions, stated as 

lim V (L) — { ^ insulator, . . 

l^oo [finite perfect metal, ^ ' 

we can see that we must go to the limit of infinite systems to distinguish 
metals from insulators. Indeed, the problem of taking the infinite-size limit is 
also occurs in the study of dielectric breakdown in Mott insulators as we shall 
see later. 



Short-time behavior — an exact diagonalization result 

Since the time evolution of many-body systems cannot be treated analyt- 
ically, we employ numerical methods to time-integrate in two steps — for 
short-time behavior and long-time behavior. For the short-time evolution 
in dielectric breakdown of Mott insulators we exactly diagonalizc the time- 
dependent Schrodinger equation as follows: First we start from the ground 
state of H((f> = 0) at time t = 0. The wave function evolves with the phase 
that increases as 

4>{t) = -> FLt/h. (45) 

Here F = eaE is the field strength, L the length of the chain. We numerically 
solve the time-dependent Schrodinger equation, 

ij t m))=Hm))\m)- m 

We choose the initial state to be the ground-state |0) of H(0), which is ob- 
tained here by the Lanczos method. The time integration of the state vector, 
which, being a many-body state, has a huge dimension, requires a reliable 
algorithm. So we adopt here the Cranck-Nicholson method that guarantees 
the unitary time evolution, where the time evolution is put into a form, 

which is unitary by definition. Here the time step is taken to be small enough 
(dt = 1.0 x 10~ 2 with the time in units of fi/t hereafter) to ensure convergence 
for L < 10, for which the dimension of the Hamiltonian is ~ 10 4 . We have 
concentrated on the total S z = subspace with — = L/2. 

Evolution of the total current 

We first plot in FigfLWb) the result for the expectation value of the current 
density averaged over the sites, J = —j; J2i a fe l T^(*) c\ +1<T a a — h.c.^j . The 
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Fig. 12. (a) The sample geometry, where an AB flux, <f>(t) = LFt, increasing 
linearly with time induces an electric force through Faraday's law. (b) Time evolution 
of the current, J(t), for a half-filled, 10-site Hubbard model for various strengths 
of the Hubbard repulsion, < U/W < 5 (W = 4t h o P is the non-interacting band 
width) for a fixed electric field F = 1/10L. Time is measured in units of Tt = h/tho P , 
LF in thop, and J(t) in l/r t . The range of the time in this panel is restricted to a 
range of the AB-flux < <f> < 1. (c) A wider plot of the current for various values of 
F with a fixed U/W = 0.25, again for the half- filled case. Here the horizontal axis 
is 4>. (d) A plot similar to (c) for a non-half- filled case (Af = Aj = 3 < L/2 — 5). 



behavior of J(t) for a fixed value of the electric field F is seen to fall upon 
three regimes when U is varied: A perfect metallic behavior ( J(t) oc t) when 
the electrons are free (U/W = 0), an insulating behavior (J(t) = 0) when 
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the interaction is strong enough (U/W 3> 1), and an intermediate regime of 
U/W where J is finite with some oscillations for finite systems. In contrast, 
a non-half-filled system in nonequilibrium (F ^ 0) has a time evolution that 
is distinct from the ground-state behavior (Fig ll'if d)). The difference has its 
root in the spectral property as will be discussed later. 

If we look at the behavior over several periods (0 < (j> < 10) for a fixed value 
of U/W for the half-filled (FigCjjJc)) and for a non-half-filled case (FigCCJd)), 
the result may be summarized as follows: 

(i) Small F regime (Mott insulator preserved at half filling) 

A drastic difference between the half-filled and doped systems appears for 
small F. When half-filled, J(t) in the limit of F — > smoothly approaches 
a periodic saw-tooth behavior with periodicity <j> = 1, which is the AB- 
oscillation of the ground-state current. 

(ii) Moderate F regime (metal) 

In this regime, the current in the half-filled case is non-zero and shows 
oscillatory behaviors (seen typically in data for LF = 0.008 in FigfTWc)). 

(iii) Large F regime (perfect metal) 

When the electric field F becomes large enough, the system behaves as 
a kind of metal. The current J(t) exhibits a long-period (A<& — <PqL) 
oscillation, which is the Bloch oscillation, a hallmark of a metal. 

The averaged current, 

{j)= ^L {j{t))dt ' (48) 

integrated over a quarter of the Bloch period (with 4>{T) = h) is plotted 
against F in Fig[l3] for various values of U. We can see that (J) becomes 
nonzero rather abruptly at the metallization as F is increased, where the 
threshold electric field increases and the F- dependence becomes weaker when 
we increase U/t. Just after the metallization some oscillation (in the F- de- 
pendence this time) is seen for finite systems. 

Evolution of the survival probability 

In order to calculate the decay rate introduced above, we compute the tem- 
poral evolution of the ground-state survival probability, 

P (t) = \(0;4>(s)\fe-*fo ffW(s))ds |0;0)| 2 , (49) 

where \0',<j>) denotes the ground-state of H (</>). The survival probability is 
related to the decay rate of the ground-state by Po(t) = e~ rt . 

The short-time feature in the survival probability is expected to be de- 
scribed by the single Landau-Zener transition between the ground-state and 
the lowest excited state, displayed by the two bold lines in the figur^l, that 

2 In Fig |14f aL three states appear in the circle. However, transition from the 
ground-state to the middle state is forbidden by symmetry. 
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Fig. 13. Dependence of the averaged current (J(t)) on F for various values of U/t 
for the half-filled Hubbard model with L = 6. 
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Fig. 14. (a) Spectrum of the half-filled Hubbard model H{<j>) for < <j> < 1. Bold 
lines represent the ground state and the first state among the tunneling-allowed 
excited states, respectively, (b) U /W-dependence of the Mott gap AE, encircled in 
(a). W = 4tho P is the non-interacting band width, and system size is L — 10 and 
U/t hop = 1.5. 



takes place around <p = I (Fig fT4l a)). If we concentrate on the two levels, 
the time evolution operator at t = At {At = -Jr is denned as the time when 
4>{At) = 1 is reached) is approximated by a 2 x 2 matrix, 

^ = ^)=( Vr ^| e " X V^fe^)' (50) 
where the tunneling probability p is given by the Landau-Zener formula [381 

ESSIES], 

r -^,U*£.\ *-JMQ/?t (5i) 



Here, A C (U) is the excitation gap (Fie l 14( b)). v = 2thop, and x the sum of 
dynamical and Stokes phases. 
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In order to verify eq. (|5ip . we can numerically calculate the survival prob- 
ability P (t) from t = to t = At for various U and F (FigfTS]). After 
determining the tunneling probability from p = 1— Pg^At), we plot it against 

the diabaticity parameter 52 • The data points (FigQIJb)) for various 
values of U fall around a common line, which is just the prediction of the 
Landau-Zener formula. The agreement is better for smaller values of U where 
we can treat the Umklapp term as a perturbation. 

Long-time behavior — a time-dependent DMRG result 

The conclusion obtained in the previous section with the exact diagonaliza- 
tion is that the short-time behavior after the electric field is switched on is 
dominated by the single Landau-Zener transition between the ground state 
and the first excited state. However, several important questions remain, e.g., 

Will the first transition remain finite in the infinite-size limit? Indeed, Kohn's 
criterion (|44f asserts that the 4> dependence of the ground-state energy of 
a Mott insulator should vanish for L — > oo. This implies that the adia- 
batic flow (Fig. fTTT c)) should become flat in this limit, which may seem to 
indicate that the transition will be washed out. However, this contradicts 
with the expression for the threshold -Fzener = ^ Ac ^l 2 ^ (eq. lfoTj) ). which 
remains finite in the L — > oo limit. Since this expression is obtained in 
a small system and in the small U limit, there is a possibility that this 
breaks down. Surprisingly, we shall show that this expression survives in 
large systems even when U is not small (Fig. [T8|). 

The effect of pair annihilation After the first transition, we expect the system 
to undergo further transitions to higher-energy levels. This process, how- 
ever, should be couterbalanced by another processes, the pair annihilation 
of doublons and holes. These processes, which do not conserve the total 
momentum in general, are caused by the Umklapp scattering. Thus the 
pair creation (= Landau-Zener transitions to high-energy states) tends to 
be offset by pair annihilation, which implies that the decay rate of the 
ground state may become smaller compared to the single Landau-Zener 
transition cas^f|. 

These questions have motivated us to study the dielectric breakdown in the 
half-filled Hubbard model for longer time periods, which is accomplished by 
the time-dependent density matrix renormalization group method. A version 
of the real-time DMRG was first intruduced by Cazalilla and Marston with a 
truncated DMRG Hilbert space and a renormalized Hamiltonian [T^]. Preci- 
sion of their method degrades rapidly in the long-time limit, since an update 
of the Hilbert space is lacking. Recently, Vidal proposed an improved method 
for simulating time-dependent phenomena in one-dimensional lattice systems 

3 A similar problem has been studied from a general point of view by Wilkinson 
and Morgan [64| . 
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Fig. 15. (a) Short-time temporal evolution of the survival probability Po(i) in the 
half-filled Hubbard model (L = 10, Nf = N± = 5) for various values of F with 
U/t = 0.1. The inset shows the solutions of the LZS equation with its asymptotic 
values indicated as dashed horizontal lines. (b)The transition probability p, with 
P(t — At) = 1 — p, plotted against the diabaticity parameter. The dashed line is 
the prediction of the Landau-Zener formula. 



employing the Trotter-Suzuki decomposition [|)TJ [55] • White and Feiguin [fj5] 
as well as other groups [18 modified this idea and combined it with the finite- 
size DMRG algorithm. 

If we denote the DMRG wave function as 
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|#>= E ^ aj+ir \l}\a 3 )\a 3+1 )\r), (52) 

where \r) is the basis of the truncated Hilbert space with dimension m 
and \aj), |<x,+i) arc the two sites that bridge the left and right blocks in the 
DMRG procedure. By employing the Trotter-Suzuki decomposition, 

-idtH -idtH 1 /2-idtH 2 /2 -idtH 2 /2 -idtH 1 /2 /ro\ 

6 — 6 6 ... 6 6 , t^OO J 

we can apply the time-evolution operator e - ldtH i/ 2 to the j-th wave function 
as 

(V^/V) = V (e-^/A i> W . a >r. (54) 

After applying e - ldtH j/ 2 ^ we diagonalize the density matrix and move to the 
next link just as in the usual finite-size algorithm. One cycle of this procedure 
results in an evolution of time by dt, and we can repeat it as many times 
as we wish. Compared with the version by Cazalilla-Marston [12] . this algo- 
rism has higher precision and we can simulate non-equilibrium excited states 
efficiently [1 8 , although one drawback of the t-dependent DMRG is that we 
can only treat systems with open boundary conditions. 

Here we study transient behaviors of the one-dimensional Hubbard model 
with open boundary condition. We use the time-independent gauge, for which 
the Hamiltonian is 

H{F) = -t hop J2 (4+i- c ^ + h - c -) +uY f n j7 n jl +F£, (55) 

i,cr j 

where X = J^j 3 n j 1S the position operator representing the tilted potential. 
As in the previous section, we start the time evolution from the F — ground- 
state |0) obtained by the usual finite-size DMRG. The wave function in this 
gauge is simply 

\&(t)) =e~ UH ^\0), (56) 
which is obtained with the t-dependent DMRG. 

Evolution of the charge density 

We first discuss the temporal evolution of the charge density, rijif) = 
(<j?(t)\nj\<I r (t)), after the electric field is switched on at t — 0. At half-filling 
the initial distribution is rij (t) = 1. After the application of the electric field, 
a charge density wave (CDW) pattern is formed when the electric field is not 
too strong (Fig. riBT a)). This state is stationary and the density profile do not 
change any further. The pattern is formed because the boundary condition 
breaks the translational symmetry, where the amplitude of the pattern corre- 
sponds to the polarization AP(F) induced by the field. When the electric field 
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Fig. 16. Temporal evolution of the charge density n, (t) in the half-filled Hubbard 
model with m = 150, L = 50, dt = 0.02 for F = 0.1(a) and F = 1.0(b). (c) depicts a 
cross section of (a) for t = 4. 

becomes stronger, charge transfers start to occur, with charge accumulation 
and charge depletion being formed around the edges in an open-boundary 
chain (Fig fl6|) . This is a sign that the ground state collapses due to quantum 
tunneling. 

The decay rate of the ground state 

The groundstate-to-groundstate transition amplitude is, in the time-independent 
gauge, 

~(t) = (0\e-* T ( H+F *)\0)ex tEo , (57) 

where we denote the ground state of H as |0) and its energy as Eq. FigurefTTTa) 
shows the temporal evolution of the ground-state survival probability |i: (i)| 2 
for a system with C//ihop = 3.5. As time evolves, the slope of — In ^(t)! 2 (oc 
the decay rate) decreases after an initial stage, which implies a suppression of 
the tunneling from the short-time behavior. This should indicate that charge 
excitations are initially produced due to the Landau-Zener tunneling from the 
ground state to the first excited states, but that scattering among the excited 
states become important as the population of the excitations grows. In other 
words, pair annihilation of carriers becomes important and acts to suppress 
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Fig. 17. (a) Temporal evolution of the ground-state survival probability |S(t)| 2 
after the electric field F is switched on at t — in the ID half-filled Hubbard model 
with U/thop — 3.5, obtained with the time-dependent DMRG for L — 50 with the 
size of the DMRG Hilbert space m = 150 and the time step dt — 0.02. The dashed 
line represents — ln|,S(t)| 2 = r(F)t + c for F/t hop = 0.17, while the dotted line 
delineates the initial slope (the short-time behavior), (b) The decay rate versus F 
in the half- filled Hubbard model. Dashed curve is a fit to eq.(l58 
is the threshold. 



where F^ ott (U) 



the tunneling rate. We have determined r(F) from the long-time behavior 
with a fitting - In \S{t)\ 2 = r(F)t + const. 

The decay rate per length r(F)/L is plotted in Fig fTTT b). where we have 
varied the system size (L = 30, 50) to check the convergence. r{F)/L is seen 
to remain vanishingly small until the field strength exceeds a threshold. To 
characterize the threshold F t h{U) for the breakdown we can evoke the form 
obtained above for the one-body system. The formula (eq. (|28l) for d = 1 
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with the error function ignored and the factor of 2 recovered for the spin 
degeneracy), 



IF 

r(F)/L = -—a(U)\n 



Mott/ 

(58) 



is originally derived for one-body problem, and an obvious interest here is 
whether the formula can be applicable if we replace the one-body F^ nd with 
the many-body F^ ott (U). In the above we have added a factor a(U), a pa- 
rameter representing the suppression of the quantum tunneling. The dashed 
line in FigfTTTb) is the fitting to the formula for £7/ihop = 3.5, where we can 
see that the fitting, including the essentially singular form in F, is surprisingly 
good, given a small number of fitting parameters. The value of a(U) turns out 
to be close to but smaller than unity (taking between 0.77 to 0.55 as U/t is 
increased from 2.5 to 5.0). 

If we perform this for various values of U we can construct a "noncqui- 
librium (dielectric-breakdown) phase diagram" , as displayed in Fig |181 which 
plots the U dependence of F^ ott . The dashed line is the prediction of the 
Landau-Zener formula |58J, 

For the size of the Mott (charge) gap we use the Bethe-ansatz result. [16] 



U J x sinh(27rytho P /t7) 

with v/th op — 2. As can be seen, the DMRG result and the Landau-Zener 
result agrees surprisingly well. 



2.5 Long-time behavior and a mapping to a quantum random walk 

Since many levels should be involved in the above pair creation/annihilation 
processes, next thing we want to have is a statistical mechanical setup for the 
time evolution of the Mott insulator. The problem at hand is a closed quantum 
system in external driving forces (e.g., electric fields), which are represented 
by a time varying parameter <f>(t) of the Hamiltonian. We want to discuss the 
asymptotic solution of the time-dependent Schrodinger equation 

inj t \*{t))=Hm))\n))- (si) 

We introduce |n; </>) as the set of eigenstates of the time-dependent Hamil- 
tonian H(4>), and denote the energy eigenvalue as E n (<f>), i.e., H (<f))\n; (/>) = 
E n ((f))\n; 4>) (Fig[T9|). Since \n;cf>) forms a complete orthonormal basis, the 
wave function ^(t)) can be expanded as 
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Fig. 18. The dielectric-breakdown phase diagram on the (U, F) plane for the one- 
dimensional Hubbard model. The symbols are the threshold F^ ott (U) obtained by 
fitting the decay rate r(F)/L to eq.([58}), while the dashed line is the prediction 
F = Fa?(U) of the Landau-Zener formula eq.j59}. 



|#(t)) - M n > ^ J ° Enis)ds \n; m (62) 

n 

with coefficients ip(n,t) = {n-^{t)\fe^^ HWs))ds \Q-Q)/e' i ^ En(s)ds . Note 
that we have removed the contribution from the dynamical phase J Q E n (s)ds/h 
in the definition of ip(n,t). Although the evolution depends on the detail of 
the system (H(<f>)), we can deduce some universal features that depend only 
on the feature of the energy levels, i.e., distribution of level repulsion in the 
spectrum. 

Each energy level is subject to the Landau-Zener tunneling to neighboring 
levels in a time period At/2, and is most conveniently expressed in terms of 
the transfer matrix representation[7l 143] . To this end, we denote the pairs as 

*(„,r) = ( H n ' T l), (63) 



and the time evolution "rule" can be expressed as 

!P(n, r + 1) = P n+1 <P(n + 1, r) + Q n -^(n - 1, r), (64) 
where P n (Q n ) is the upper (lower) half of a 2 x 2 unitary matrix, 

°>-(zty «■=("» »)■ "--(ID- <65) 

The diagonal elements of U n represent the Landau-Zener transition from 
the n-th level to (n — 1) or (n + l)-th level, where the explicit form is 



30 Takashi Oka, Hideo Aoki 



(a) Half-filled (Mott insulator) n (b) idealized spectrum 




Fig. 19. (a) The spectrum of the half-filled Hubbard model. The circle corresponds 
to Landau-Zener transition between the two energy levels which can be expressed by 
a 2 x 2 unitary matrix (eq. (|66p ). (b) Idealized energy levels where level anti-crossings 
are expressed by circles. Quantum interference takes place when contributions from 
different paths are considered. 




t t+1 



Fig. 20. Application of Q n and P„ in the quantum random walk. 



Here the Landau-Zener tunneling probability p n depends on the ratio of the 
Zener threshold field i*z ener and the electric field F as 

Pn = CX P ( -TT^Lj , ( 67 ) 

where F% cnGT generically depends on n. 

^ ) n ' A J as a "qubit" on "site" n, ea. fM]) defines 

an evolution of a one- dimensional quantum walk with a reflecting boundary 
at n = corresponding to the ground-state fFig rTW cM. A quantum walk is 
a quantum counterpart of the classical random walk. Models with essentially 
equivalent ideas have appeared in various fields: to name a few, quantum 
transport and dissipation [TTJ [5H] , quantum Hall effect [13] , optics [S3 [S] and 
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recently in quantum information @51 [331 IM 1121 IM 1111 EOl HOI [Ml HDJ [B] . 
In the field of quantum information (see e.g. [M]), introduced by Aharonov, 
Ambainis, Kempe and Vazirani in 2001 [57], the quantum walk is arousing 
interest in hope of revealing new features in the quantum algorithms (for 
reviews see [28] [40l [6] ) . Researches stem into many directions, e.g., the effect 
of absorbing boundary conditions [22l [30l [53], higher-dimensional systems 
[54] [40] [32] , localization in systems with internal degrees of freedom [21] , and 
many powerful analytical techniques are being developed. 

An important feature of the quantum walk, as opposed to the classical 
walk, is that different transition paths interfere with each other quantum 
mechanically. We in fact find that the quantum interference leads to a dy- 
namical localization, an analog of Anderson's localization taking place in the 
energy space rather than in the position space. In our previous work [Jj we 
have employed the PRQS method, a technique to treat quantum walks, to per- 
form the path integral, and obtained the exact asymptotic distributions of the 
wave function for a simplified model. The resultant states can be categorized 
in three types depending on the strength of the electric field, as schemati- 
cally plotted in Fig. [21] (A) is an adiabatic evolution that takes place in weak 
driving forces (electric fields). The dominant part of the distribution p(n,t) 
against energy is a delta function localized around the ground-state. When the 
driving force become stronger, quantum tunneling broadens the delta func- 
tion, as plotted in (B). The shape of the peak is maintained by a balance 
between tunneling and dynamical localization. (C) is the case where the driv- 
ing force overwhelms the effect of dynamical localization, and the system is 
driven rapidly into the excited states. This, in our view, corresponds to the 
dielectric breakdown. 
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Fig. 21. The time evolution of the distributions of wave function amplitude 
p{n,i) = \ip(n, t)\ 2 in energy space. The vertical axis n is the index of the energy 
levels. (a)For small tunneling, the distribution is localized at the ground state, (b) 
For intermediate tunneling, a localized state remains, while the amplitude starts to 
bifircate into excited states of the wave function is excited, (c) When the tunneling 
is larger than the threshold, the localized state disappears. 



2.6 Experimental implications 

Now we discuss experimental implication of the many-body Landau-Zener 
transition mechanism. In fact, there are several mechanisms which may lead to 
breakdown of insulators. For example, Frohlich's electric avalanche mechanism 
may take place, in which a small number of excited electrons act as a seed 
and become accelerated by the electric field until they cooperatively destroy 
the insulator. We can distinguish Landau-Zener transition from the avalanche 
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mechanism through the temperature dependence and from interface effects 
by changing the size of the sample. Another important effect is the band 
bending near an interface of the Mott insulator and electrodes. For a thin 
sample, this may lead to injection of carriers, and results in the interface 
Mott transition! 591. 



(a) 



oscilloscope 



Vaivage 




Source-measure unit 




E (V/cm) 



Ohmic law 
(finite T effect) 

Fig. 22. (a) A schematic measurement circuit, (b) J— E characteristics of a SrCu02 
sample, (c) Temperature dependence of the resistance for various applied voltages. 
After Taguchi et al. gT]. 
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Dielectric breakdown of one-dimensional Mott insulators was experimen- 
tally studied by Taguchi et ai, who obtained the J—E characteristics (Fig l^)) 
of Sr 2 Cu03 and SrCuC>2 samples (31). which are both quasi-lD, strongly corre- 
lated electron systems. Experiments were done by placing small single crystals 
in circuits as shown in Fig|22Fa). and the voltage drop V was measured while 
the current density J was fixed. Depending on the strength of the electric field 
F, transport properties change drastically, as summarized in the following. 

In weak electric fields, the J—E characteristics shows an Ohmic behavior 
at finite temperatures. When the electric field exceeds a threshold value, the 
current shows a dramatic increase. Such drastic changes cannot be explained 
by perturbation in F, and we must consider non-perturbative effects, i.e., a 
behavior essentially singular in F like J ~ function ofexp (—F t h/F) , which 
is a typical tunneling effect with threshold i^h- The temperature dependence 
(FigllSD of the threshold can be fit well by F th (T) / F th (0) ~ exp (-T/T ). This 
excludes the avalanche mechanism, for which an activation type temperature 
dependence (i^ alanche (T)/.F th (0) ~ exp(T /T)) is expected. 

One indication that the breakdown is indeed quantum in nature is that 
the threshold extrapolates to a finite value for T — ► 0. From the extrapolation 
(FigUHfa)), we obtain a threshold, 

F t c h xp ~ 10 6 - I0 7 (eV/cm), (68) 

for Sr2Cu03 and SrCu02- The Landau-Zener result (intended for T = 0) of 
the threshold (eq. J5TJ)) is 

\A (77V21 2 

(P) = ~ (leV) 2 /(10- 7 eV/cm) ~ 10 6 (eV/cm) (69) 

is comparable with the experimental result. 

Interestingly, the decay rate r(U) we have introduced theoretically can 
be measured experimentally 41 . This is done by studying the transient be- 
havior of the current after the electric field is switched on at t = 0. At first 
the current density is zero, and then becomes non-zero after a certain delay 
time t = t( F)(FigfI3[b), lower inset). The authors in [IT] have introduced a 
phcnomenological percolation model to relate the delay time with the produc- 
tion rate P(F) of the conductive domains (see [H]). In this model, conductive 
domains are envisaged to grow in the sample, and the current density is as- 
sumed to become finite when the left and right electrodes are connected by 
these domains (Fig |25f c)). This leads to a relation, 

PW = -(F±y\ (70) 

The experimental result for the production rate in Fig l23f b) is obtained in 
this way. 

The nature and the microscopic origin of the "conductive domain" are not 
clear, but if we interpret them to be domains with a high density of charge 
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Fig. 23. An experimental result for the temperature dependence of the threshold 
electric field F t h{T) for the dielectric breakdown (adopted from Taguchi et al. [41] ). 
Dashed lines correspond to a fitting to exp(Tb/T) predicted by Frohlich's electron 
avalanche. Solid lines correspond to a fitting to exp(— T/To). (b) The electric-field 
dependence of the delay time Td for Sr2Cu03. (c) The production rate of "conductive 
domains". The inset exemplifies the temporal evolution of the current at various 
applied voltages. Adopted from Taguchi et al. [41] . 



excitations produced by the Landau-Zener transition, the vacuum decay rate 
per volume r(F)/L d characterizes quantum tunneling from the ground-state 
to excited states. With the identification we expect that the decay rate and 
the production rate are identical, i.e., 

P(F) ~ r(F). (71) 

This identification is encouraged by the field dependence of P(F) (lower panel 

of Figl237b)). which is close to the expected form (FjL ~ —^-a(U) In 1 - exp ( -n th F j ) 

of the decay rate. 

In this experiment a scaling study — a systematic change of the size of 
the sample — was also performed to confirm that the nonlinear effect occurs 
in the bulk. From these observations, we conclude that the experiment by 
Taguchi et al. |41j can be explained by the many-body Landau-Zener tunneling 
mechanism. However, to be more confident, we need to know the temperature 
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dependence of the threshold theoretically, which is still a challenging task in 
the present many-body system. 

2.7 Conclusion 

In this article, we have explained how dielectric breakdown of Mott insulators 
can be explained from the nonequilibrium behaviors of charge carriers, espe- 
cially from their creation and annihilation processes. Both processes are the 
result of many-body Landau-Zener nonadiabatic tunneling transition between 
many-body energy levels, where charge creation processes are counterbalanced 
by annihilation processes. From numerical result we have obtained a nonequi- 
librium (dielectric-breakdown) phase diagram. If the coherence of the dynam- 
ics is preserved at sufficiently low temperatures, a quantum interference, as 
modeled by a quantum walk in energy space, may lead to dynamical local- 
ization, which saturates the creation process and leads to a non-equilibrium 
stable state. The decay rate r(F) that we have discussed is a measurable quan- 
tity: it is the production rate observed by Taguchi et al. [41j in copper oxides. 
The experimental result is consistent with our prediction r(F) ~ £ ln(l -p) 
when an extrapolation to zero temperature is made. It is an interesting future 
topic to understand the properties of the non-equilibrium stable state in more 
detail. 

An important open question is how the energy dissipation processes take 
place in nonequilibrium situations. Here we have stressed that the many-body 
processes act effectively as a source of dissipation through scattering, but an 
explicit incorporation of heat-bath effects, electrode effects, etc, is left to a 
future problem. 
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